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ABSTRACT 

We find that the majority of systems hosting multiple tidal disruptions are likely to contain hard 
binary SMBH systems, and also show that the rates of these repeated events are high enough to 
be detected by LSST over its lifetime. Therefore, these multiple tidal disruption events provide a 
novel method to identify super-massive black hole (SMBH) binary systems with parsec to sub-parsec 
separations. The rates of tidal disruptions are investigated using simulations of non-interacting stars 
initially orbiting a primary SMBH and the potential of the model stellar cusp. The stars are then 
evolved forward in time and perturbed by a secondary SMBH inspiraling from the edge of the cusp to 
its stalling radius. We find with conservative magnitude estimates that the next generation transient 
survey LSST should detect multiple tidal disruptions in approximately 3 galaxies over 5 years of 
observation, though less conservative estimates could increase this rate by an order of magnitude. 
Subject headings: black hole physics — galaxies: kinematics and dynamics — galaxies: evolution — 
galaxies: nuclei 



1. INTRODUCTION 
Stars with radius rv, and mass M*, which pass within 
the tidal disruption radius r t ~ ^(Mbh/^*) 1 ' 3 01 a 
super- massive black hole (SMBH) of mass Mbh, will be 
ripped apart by tidal forces. In the case of sun-like stars, 

r t « l.2r s M~ 2/3 , (1) 

where r s is the Schwarzschild radius, and Ms is 
M BH /1O 8 M . Therefore, when M BH £ 10 8 M o the 
Schwarzschild radius lies outside r t and any sun-sized 
star would be swallowed whole. Below this critical black 
hole mass the star's debris is launched on orbits whic h 
span an energy range AE w GM B w\/r 2 t (iReesl H988I) . 
This energy range is large compared to the energy of 
the highly elliptic initial orbit, and hence half the ma- 
terial will be unbound while half will fall back onto the 
black hole. For main sequence stars the fall back rate 
declines as i _5/3 ()Phinnevlll989| ). This fall back rate is 
initially su per-Eddington for the cano nical 10% accretion 
efficiency ( (Evans fe Kochanekl Il989f ). but it is unclear 
whether the disk adjusts to lower its accretion rate or 
a rad i atively driven outflow results (|Strubbe fc Quataertl 
l2009HLodato fc Rossill2010h . 

Galaxies harboring an isolated SMBH at their center 
are expected to quickly clear a 'loss cone' of orbits whose 
angular momenta about the black hole are low enough 
that their peribothra lie inside r t . At this point tidal dis- 
ruptions are predicted at a rate ~ 10~ 4 — 10~ 5 yr~ 1 as 
stars diffuse into the loss cone dMagorrian fe Tremainej 
[19991: IWang fc Merrittl l200i IDonlev et al.l I2002TI . The 
majority of candidate tid al disruptions thus f ar have been 
found th ough X-ray (e .g. IDonlev et al.ll2002r i or UV sur- 
veys (e.g. lGezari et al.ll2008lh This is expected, as can be 
seen by modeling the tidal disruption as a thick disk emit- 
ting as a black body with luminosity L c dd, temperature 
T e ff, and initially extending to r t . In reality the disk will 
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expand outwards on a viscous timescale and the initial 
super Eddington rate could launch an outflow. I gnoring 
these complications, however, gives (lUlmerlll999 h 

/ , r \ -1/6 / S, -1/2 



T eS ~ 3.7 x 10 M 



and the spectrum peaks in the extreme UV. Despite be- 
ing in the Raylcigh- Jeans tail of this flux, optical tran- 
sient surveys such as the Palomar Transient Factory 
(PTF), the Panoramic Survey Telescope and Rapid Re- 
sponse System (Pan-STARRS) and the Large Synoptic 
Survey Telescope (LSST) provide the prospect of finding 
many more tidal disruptions because of their unprece- 
dented combination of high cadence and depth. It is ex- 
ected that LSST will detect a striking ~ 100-3000 yr^ 1 
Strubbe fc Quataertl 120091) . where the major uncer- 
tainty is how the luminosity from the uncertain super- 
Eddington phase of the tidal disruption is included. 

In this letter we calculate the rates of multiple tidal 
disruptions from the same merging SMBH binary system, 
and show that the detection of multiple tidal disruptions 
from a single galaxy likely indicates the galaxy hosts a 
SMBH binary with a parsec to subparsec separation. Our 
results are summarized in table [1] 

2. SIMULATIONS 

Two mechanisms for enhanced rates of tidal disrup- 
tions in close SMBH binary systems ha ve been con- 
sider ed in the literature, the Ko zai effect (|Iyanov et ail 
12001 ) and chaotic 3-body orbits (|Chen et al.ll2009D . We 
have extended this work by performing a series of re- 
stricted 3-body scattering experiments including two ad- 
ditional key aspects: the evolution of the binary, and 
the non-Keplarian stellar potential. Our simulations are 
described below, where we use G = 1 throughout. 

In considering the stellar potential, it is important 
to have a stellar distribution consistent with a central 
SMBH. Thus, the primary SMBH of mass Mi was placed 
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in a Tremaine-model (also known as the Dehnen-modcl) 
cusp ()Tremaine et al.l ll994: Dehncn 1993) whose density 
is given by 



p(r) 



1 



4-7T r 3 - r '(l + r) 1 + f > 



(3) 



The advantage of this model is that it self-consistently 
describes a finite mass of stars distributed around a cen- 
tral black hole together with their stellar potential 
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rj = l, 



(4) 



where p is the ratio of Mi to the total stellar mass in 
the cusp. Throughout this work we use p = 0.5. Since 
the stellar mass is normalized to unity this results in 
the primary's mass being half the cusp's stellar mass. A 
Bachcall-Wolf cluster corresponds to r\ = 1.25. 

In this way our simulations for a given rj and q then 
depend only on the tidal disruption radius rt . The scaling 
of these parameters to real galaxies is described later by 
equations [7] and [8] 

Each star's initial position was chosen according to 
equation [3] with a random orientation. To pick each 
star's velocity, we numerically tabulated the distribu- 
tion function f(£) and then picked an isotropic ve- 
locity from the required distribution, 4irf(£)v 2 dv = 
47r/(£ )y/2(^(r) — E) d£, using the rejection method. 

The stars are initially on an orbit consistent with the 
primary SMBH and the stellar potential, however their 
orbits are perturbed by the secondary SMBH whose orbit 
we evolve with time approximately following an inspiral 
dominated by dynamical friction. To model this inspiral, 
the secondary SMBH is initially given an eccentricity of 
0.1 and a binary separation equal to the cusp radius, r c . 
It is then migrated inwards on a path governed by 



~dt 



[p(l + q) + M ir (<r)] 



w df v 



(5) 



where Af*(< r) is the stellar mass interior to r, q < 1 is 
the binary mass ratio and 



47rlogA(7 p p(< v) 



(6) 



chara cterizes the dynamical friction (|Binnev fc Tremaind 
l2008f l. Here p(< v) is the density of stars at r with veloc- 
ity less than v. We have used a Coulomb logarithm that 
begins at log A sa 4, but which smoo thly decreases to zer o 
at the stalling radius calculated by iSesana et al.l (|2008t ). 
The functional form of the decrease was chosen to ap- 
proximate the rate of shrinkage caused by the ejection of 
stars during our scattering experiments. For instance, in 
the lower panels of figure [1] we plot the change in stellar 
and binary energies. If the functional form of log A had 
been chosen perfectly the two would lie on top of each 
other. 

To perform our scattering experiments we have imple- 
mented the adaptive sym plectic integrator described in 
iPreto fc Tremainel ()1999j ). with a timestcp which varies 
as At ~ U^ 1 where U is the potential energy. With 
this choice of timestep the integrator has the desirable 




Figure 1. Results for mass ratio q = 0.1. The left hand panels 
show an r\ = 1.25, Bahcall-Wolf cusp, the right hand panels an 
rj = 1.5 cusp. The upper panels shows tidal disruption rate, IY;,,, 
for tidal disruption to cusp radius ratios rt /r c = (9, 7, 5) X 10 
in solid, dotted and dashed lines respectively. The distribution of 
disruption times have been kernel smoothed with a Gaussian of 
width <j = 2. The middle panels shows the evolution of the binary 
separation as a solid line and radii enclosing 0.1, 0.2, 0.4% of the 
stellar mass in dotted, 1, 2, 4% in dashed and 10, 20, 40% in dash- 
dot lines. The lower panel shows the evolution of the energy of 
the binary as the solid line and the stars as the dotted. In a fully 
self-consistent evolution these would lie on top of each other. The 
simulations are scaled by fl = (2GMi/rf) 1 ^ 2 . 
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Figure 2. Results as described in figure [T] for q = 0.3. 

property that it reproduces exactly Keplarian orbits in 
Keplarian potentials, allowing Kozai resonances to be 
correctly reproduced since the spurious precession fre- 
quently found in other algorithms is absent. In addition, 
this integrator is well suited for this problem since it has 
been shown to correctly reproduce the highly ecc entric 
orbit s required for a star to be tidally disrupted (jPeterl 
[2001 . In our simulations stars are considered disrupted 
when they pass within the tidal disruption radius of ei- 
ther hole appropriate for a sun-like star. Wc find the vast 
majority of stars are disrupted by the primary SMBH. 

It should be noted that for the duration of the scat- 
tering experiments the stellar potential is assumed to be 
centered on the primary SMBH and is not allowed to vary 
with time. These assumptions are made for simplicity, 
and will be relaxed in future studies. 
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Figure 3. Plots of our rj = 1.25, q = 0.3 simulation. Left panel 
shows the stars that are tidally disrupted for rt/r c = 5 X 10 — 7 
as a function of their in initial radius and z-component of angular 
momentum normalized to the circular angular momentum at that 
radius, J z /J c . The Kozai wedge is plotted together with the over- 
all stellar density. A large fraction of the disrupted stars lie well 
outside the Kozai wedge indicating that these are chaotic orbits (cf. 
Chen ct al. 2009). The contours show the initial stellar distribu- 
tion, each is evenly spaced in density. The right hand panel shows 
the rates from the same simulation scaled using the relations in 
equations [7] and [8] for Ah = 10 8 M Q in solid, 5 X 1O 7 M in dotted 
and 10 7 Mq dashed lines. 

The results of our scattering experiments and the re- 
sultant disruption rates, r t d, are shown in figures Q] and 
[2] for mass ratios q = 0.1 and q = 0.3, respectively. The 
plots are scaled by Q = {2GM 1 /rl) 1 / 2 . 

We find that the majority of tidal disrupti ons are due 
to the type of chaotic orbits described by iChen et all 
(120091) as opposed to the Kozai effect discussed by 



(I^UUiH) as opp oE 
llvanov et all pO 



20051) . This is demonstrated in the left 
hand panel of figurc[3J The reason for this is twofold: pre- 
cession of the orbits of the stars in the non-Kcplarian po- 
tential destroys the secular Kozai effect for the majority 
of orbits, and the few stars closely bound to the primary 
whose precession rate is lower have Kozai timescales 
longer than our simulation. 

Our rates are lower than those discussed by I Chen et al.l 
(|2009t ) largely because we have considered less steep 
cusps. This both reduces the number of stars that can 
be disrupted as the binary hardens, and increases the 
orbital timescale at the hardening radius. Both effects 
reduce the rate of disruptions. In addition we have con- 
sidered the binary evolution which IChen et all (|2009D did 
not, although this has a smaller effect 

To apply o ur simulations to phy sical galaxies, we use 
the fits from iMerritt et all (120091) to the i nner r egions 
of ACS Virgo Cluster galaxies (|C6te et al.l I2001 . For 
power-law galaxies these giveQ 

r. = 22 (A/i/lO^o) - 55 pc , (7) 

where r, is defined such that the stellar mass interior to 
r. is 2Mi and Mi is the mass of the SMBH. Matching this 
to the Tremainc model such that the central densities are 
equal gives r c = r m for fi = 0.5, used in our simulations. 
The ratio r t /r c is 

r t /r c = 4.9 x lO" 7 (Mi/lO 8 M ) a22 . (8) 

With these searings, our simulations for i] = 1.25 and 
q = 0.3 are shown in figure [3] 

2 D. Merritt person al communication. From fitting to figure 2 
of IMerritt et all 1(20091 ). 



3. OBSERVABLE TIDAL DISRUPTIONS 

The absolute magnitude of an individual disruption is 
likely to depend on complex details with considerable 
uncertainties in modeling quantities such as the SMBH 
mass, the SMBH spin and t he geometry of the disruption 
(|Strubbe fe Quataertll2009T) . 

We instead derive a simple empi rical estimate of th e 
volume accessible by comparison to iGezari et al.1 (|2008| ) . 
Two luminous optical events coincident with UV flares 
were discovered in ~ 2.9 deg 2 . Their spectra and light 
curves were consistent with tidal disruption events, mak- 
ing this their most likely explanation. Their redshifts 
were z = 0.33 and z = 0.37, giving extinction corrected 
(but not K-correc ted) absolute <?-band magnitudes of 
-17.7 and -18.9 (|Gezaril [2010h . Requiring that these 
two cases be 2 mags brighter than the 25.0 <?-band limit 
of LSST gives maximum redshifts of detection of z — 0.27 
and z = 0.43, respectively. The 2 magnitude buffer bet- 
ter ensures a convincing light curve, which would display 
the characteristic fast rise and decay of a tidal disrup- 
tion. Based on these numbers we choose z = 0.35 as the 
limit for LSST. 

There is also only a small range of SMBH masses 
which needs to be considered. Because SMBHs of mass 
greater than 10 8 M Q can't tidally disrupt stars (equation 
[J ignoring SMBH spin), and SMBHs of mass less than 
10 7 M Q are significantly less luminous, (particularly if the 
super-Eddington phase is neglected) , and so can only be 
observed in a much smaller volume, we restrict our anal- 
ysis to SMBHs with masses between 10 7 — 10 8 M Q . 

4. RATES OF SINGLE TIDAL DISRUPTIONS 

We now calculate the rate of tidal disruptions observ- 
able by LSST both for systems with isolated SMBHs and 
systems with SMBH binaries. 

In the case of isolated SMBHs, the rate of tidal disrup- 
tions observed by LSST will be 



^(td) _ „ 
'Single ~ J sky 



dN 



dNL 



-V c (M B H)rtd(MBH)dM B H (9) 



where / s k y is the fraction of the sky covered by LSST, 
Ttd is the rate of tidal disruptions per galaxy, V c is 
the total comoving volume over which a tidal disrup- 
tion is observable and dN/dMsH is the black hole mass 
function. We have used the black hole mass function 
(|Aller fc Richstondl2002l) 



dN 

dAlBH 



(MshX 



a -M BH /M* H 



(10) 



with the parameters c = 3 x 10 11 M Q 1 Mpc 3 , 
M£ H = 1.1 x 10 8 M W and a = 0.95 (values de- 
rived by lAller fc Richstond 120021 scaled to H = 
71 kms -1 Mpc -1 ). The rate of tidal disruptions per 
galaxy is highly uncertain, and so we parameterize, 
Ttd = 7 x 10 _5 yr _1 , scaling to the observationalfy mo- 
tivated constant rate per galaxy of ([Donley et al.l I2002T) 
independent of Mbh- We also assume V c is indepen- 
dent of Mbh and is 10.7 Gpc 3 , corresponding to our 
rcdshift limit of z = 0.35 with the assumption that 
Hq = 71 kms -1 Mpc -1 . Using these approximations and 
/sky ~ 0.5, the rate of tidal disruptions detected by LSST 
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in galaxies containing isolated SMBHs (equation [5]) is 
predicted to be 



•10 s Af o j N 



10 7 M Q dM B n 



BH 



:3007 yr~* 



(11) 



Now consider systems hosting binary SMBHs. The 
rate of disruptions will be 

^bin = /sky / V c T bin (M u q, t)i? merge (Mi, g) dq dM 1 dt , 

(12) 

where R meTge (Mi,q) dqdAIi is the rate of mergers per 
unit comoving volume for binary SMBHs with primary 
mass between Mi and Mi + dM\ and with mass ratio 
between q and q + dq. The quantity J Thin (Mi, q, t) dt 
is the total number of tidal disruptions in a merger and 
was linearly interpolated from the simulations in figures 
[JJ and [2] together with equations [JJ and [5] 

Over the narrow range of redshift and primary mass 
accessible we approximate 



,(M!,q) 



(13) 



where we have assumed the mass ratio distribution of 
SMBH binaries follows th e local galaxy merger mass ra- 
tio distribution given bv IStewart eTatl (pool . F(q) = 
q°- 25 (l — q) 11 . The normalization constant C was chosen 
to be 0.05 Gyr -1 , so that we we reproduce the simulated 
local merger rat<3 of SMBHs with primary SMBH mass 
between 10 7 M Q and 10 8 M© and q > 0.05 of approxi- 
mately 9 x 10~ 5 Mpc~ 3 Gyr" 1 . 

Using these approximations, then, if all galaxies have 
an i] = 1.25 cusp we can expect LSST to detect 

1O 8 M 0.5 
X J dAh J dq J dt-^-F(q)T hin (M l7 q,t) 



10 7 M o 0.05 

=10 yr- 1 , 



(14) 



where we have limited the mass ratio to q < 0.5 since our 
simulations assume that the secondary has been stripped 
of stars when it reaches the cusp. 

5. RATES OF MULTIPLE TIDAL DISRUPTIONS 

We now calculate the rate of multiple tidal disruptions 
in systems containing isolated SMBHs. Over a period of 
observing i b s the total number of tidal disruptions will 
follow a Poisson distribution with mean iobsLtd- The 
probability of observing multiple tidal disruptions from 
a single galaxy is therefore 1 — Q(2, i b s rtd), where Q is 
the incomplete gamma function. 

Then the expected number of isolated SMBHs exhibit- 

3 M. Volonteri person al communication. From data in figure 2 
of lVolonteri et al.l 1120091') . 



ing multiple tidal disruptions during i b s is 

1O 8 M 



yy(multi) _ 
single 



sky 



dN 

dAh 



BH 



[1-Q(2, to^td)] gJMbh 



10 7 A4q 



0.03 7 2 (W5yr)' 



(15) 



w (multi) 
JV bin 



Similarly, the expected number of multiple tidal dis- 
ruptions observed from binary SMBHs is 

/sky J V c [1 - Q(2, r bin t obs )] fl m erge dAh dt dq . 

(16) 

Using the same approximations used in estimating equa- 
tion [14] we find over an observation time, i bs = 5 yr, 
the expected number of close binary SMBHs exhibiting 
multiple tidal disruptions observable by LSST to be 

r 10 8 M Q r 0.5 

KZ M) ~VJ sky C / dMi / dq 

J10 7 M Q JO. 05 

dN 

I dt^r^F(q) [1 - Q(2, r bin (Mi, q, t)t ohs )} 



dAL 



BH 



=3. 



(17) 



where the quantity in square brackets was calculated 
from our simulations together with equations [7] and [5J 

Towards the upper end of the range 10 7 — 10 8 M Q 
Tbiniobs ~ 0.5 for major mergers. This indicates that 
the majority of close SMBH binaries with primaries in 
the upper end of this range could potentially be iden- 
tified using multiple disruptions. The expression above 
broadly scales as t^ hs but because Tbmiobs ~ 0.5 for some 
systems this is only approximate. 

6. DISCUSSION 

We have estimated the enhanced rate of tidal disrup- 
tions from SMBH binaries, and shown that if a system 
exhibiting multiple tidal disruptions is observed then in 
our fiducial model it is <~ 100 times more likely to be a 
close SMBH binary than an isolated SMBH system. It 
has also been shown that the upcoming transient survey 
LSST is likely to detect several systems with multiple 
disruptions during a 5 year observation period. 

Once a double tidal disruption is detected these galax- 
ies would be expected to have a steady tidal disruption 
rate, with further events on a human timescale. In addi- 
tion, a double disruption makes the system immediately 
worthy of detailed study and monitoring to estimate the 
properties of the predicted binary black hole. 

We have also shown that in our fiducial model ap- 
proximately 3% of all tidal disruptions occur in binaries. 
Other signatures may identify tidal disruptions that oc- 
curred in binaries. These include possible spectroscopic 
signatures, morphology or kinematics indicating a recent 
major merger, or interruption of the tidal disru ption flare 
on a binary orbital timescale (jLiu et al.ll2009| ). 

These results depend strongly on several quantities. 
In particular, the tidal disruption rates are largely deter- 
mined by the number of stars in the central regions of 
the galaxy, which, in turn, depends on the cusp profile 
and the size of the cusp. In this sense, multiple tidal 
disruptions are also diagnostic of cusp profiles. We have 
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21.0 (g-band) 


0.06 


0.2 


0.77 


0.02 


0.9 x 10" 4 7 i! 


0.007 
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0.35 
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10 


0.03 7 2 
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Table 1 

Summary of rates of tidal disruptions for three current and upcoming transient surveys. Symbols and calculation are described in the 
text. The numbers A r '. mu ! t1 ' and 7v' mult1 ' are for an observation time of t nh „ = 5 yr and scale roughly as i 2 , . 

single Dm °us j o j obs 



assumed Bahcall-Wolf type cusps for our rates. More- 
over, the rates of multiple tidal disruptions from isolated 
SMBHs depend on the square of the uncertain tidal dis- 
ruption rate, for which we have adopted a conservative 
fiducial rate of 10~ 5 yr _1 . 

All our numbers scale by the uncertain detection vol- 
ume, which could b e significantly higher th an we have 
assumed. Recently Ivan Velzen et aIT(|2010f ) found two 
candidate disruptions with absolute g-band magnitudes 
—20.3 and —18.3. If the event with magnitude —20.3 
was representative of the higher black hole mass disrup- 
tions where our binary induced disruptions typically oc- 
cur, then LSST could detect disruptions of this type to 
z ~ 0.7 increasing our predicted rates by approximately 
an order of magnitude. 
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